/*
 * Operations on binary float-encoded symmetric matrices.
 *
 * Usage:
 * n <lhs> add <rhs>            -- set out_ij = lhs_ij + rhs_ij
 * n <lhs> mul <rhs>            -- set out_ij = lhs_ij * rhs_ij
 * n <mat> scale <t>            -- set out_ij = t * m_ij
 * n <mat> shrink <t>           -- set out_ij = m_ij / (m_ij + t)
 * n <a> shrink_with <b> by <t> -- set out_ij = a_ij*b_ij / (b_ij + t)
 * where in all cases n is the size of the matrix (number of rows, which is
 * the same as the number of columns because the input matrices are assumed
 * to be square). In the scale and shrink cases, you can pass "-" for the
 * <mat> argument, to use standard in.
 */
#include <cstdlib>
#include <fstream>
#include <iostream>
#include "u1.h"
#include "symmetric_matrix.h"

using namespace std;

int main(int argc, char **argv) {
  if (argc != 5 && argc != 7) {
    cerr<<"See comments for usage."<<endl;
    exit(EXIT_FAILURE);
  }
  size_t n = atoi(argv[1]);
  string lhs_string = argv[2];
  string op_string = argv[3];
  string rhs_string = argv[4];

  if (op_string == "add" || op_string == "mul") {
    if (argc != 5) {
      cerr<<"See comments for usage."<<endl;
      exit(EXIT_FAILURE);
    }

    symmetric_matrix<float> lhs(n), rhs(n);

    ifstream lhs_is(lhs_string.c_str());
    lhs.read_binary(lhs_is);
    lhs_is.close();
    cerr<<"Read LHS."<<endl;

    ifstream rhs_is(rhs_string.c_str());
    rhs.read_binary(rhs_is);
    rhs_is.close();
    cerr<<"Read RHS."<<endl;

    if (op_string == "add") {
      lhs += rhs;
    } else if (op_string == "mul") {
      lhs *= rhs;
    }

    lhs.write_binary(cout);
  } else if (op_string == "scale" || op_string == "shrink") {
    if (argc != 5) {
      cerr<<"See comments for usage."<<endl;
      exit(EXIT_FAILURE);
    }

    symmetric_matrix<float> mat(n);

    if (lhs_string == "-") {
      cerr<<"reading mat from stdin"<<endl;
      mat.read_binary(cin);
    } else {
      cerr<<"reading mat from "<<lhs_string<<endl;
      ifstream mat_is(lhs_string.c_str());
      mat.read_binary(mat_is);
      mat_is.close();
    }
    cerr<<"Read matrix."<<endl;

    float t = atof(rhs_string.c_str());

    if (op_string == "scale") {
      mat *= t;
    } else if (op_string == "shrink") {
      symmetric_matrix<float> den(n);
      den += mat;
      den += t;
      mat /= den;
    }

    mat.write_binary(cout);
  } else if (op_string == "shrink_with") {
    if (argc != 7) {
      cerr<<"See comments for usage."<<endl;
      exit(EXIT_FAILURE);
    }
    
    float t = atof(argv[6]);
    symmetric_matrix<float> a(n), b(n);

    ifstream lhs_is(lhs_string.c_str());
    a.read_binary(lhs_is);
    lhs_is.close();
    cerr<<"Read LHS."<<endl;

    ifstream rhs_is(rhs_string.c_str());
    b.read_binary(rhs_is);
    rhs_is.close();
    cerr<<"Read RHS."<<endl;

    // Don't use temporary storage for this. The result is: a.*b ./ (b + t)
    a *= b;
    b += t;
    a /= b;

    a.write_binary(cout);
  } else {
    cerr<<"Unknown operation "<<op_string<<"."<<endl;
    exit(EXIT_FAILURE);
  }

  return EXIT_SUCCESS;
}
/*
* Copyright (c) 2009 John Lees-Miller
* 
* Permission is hereby granted, free of charge, to any person
* obtaining a copy of this software and associated documentation
* files (the "Software"), to deal in the Software without
* restriction, including without limitation the rights to use,
* copy, modify, merge, publish, distribute, sublicense, and/or sell
* copies of the Software, and to permit persons to whom the
* Software is furnished to do so, subject to the following
* conditions:
* 
* The above copyright notice and this permission notice shall be
* included in all copies or substantial portions of the Software.
* 
* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
* EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
* OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
* NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
* HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
* WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
* FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
* OTHER DEALINGS IN THE SOFTWARE.
*/

